subroutine srccsd_io(ie, fij, fab, fia, & 
          iijkl, iijka, iabij, &
          t1, t2, r1, r2, &
          fae, fmi, fme, &
          wmnij, icore, & 
          maxcore, escf, eccsd, rivv, rivo, rioo) 
use mod_iop
use mod_size
use mod_orbit
use mod_ioff
use mod_ri
implicit none
include 'omp_lib.h'
!
integer prec2, isd2, igpu2, totgpu2
common/precsoc/prec2, isd2, igpu2, totgpu2
!
integer, intent(in) :: maxcore 
real*8, intent(inout) :: escf, eccsd
real*4, intent(inout) :: ie(*), fij(*), fab(*), fia(*) 
real*4, intent(inout) :: t1(*), t2(*), r1(*), r2(*) 
real*4, intent(inout) :: iijka(*), iabij(*), iijkl(*)
real*4, intent(inout) :: fae(*), fmi(*), fme(*)
real*4, intent(inout) :: wmnij(*)
real*4, intent(inout) :: icore(*)
real*4, intent(inout) :: rivv(*), rivo(*), rioo(*)
!
integer nsize, i, j, a, b, icrsizg
integer maxvar, maxdim, iter
integer i0, i1, i2, i3, i4, i5, i6, i7, i8, i9
integer i10, i11, i12
real*8 t1diag
real*8 dtmp, errdiis, eccsd2, ecc_i(200), Eccsdio(1000)
real*8 eso, tt(20), rt(20), tt1, tt2, iter_time(200)
real*4, allocatable :: idoso(:), itmp(:)
!
real*8, external :: ddot, dnrm2, nrm2so
integer, external :: idsymsoc, irpdso
!
ldebug = .false.
ldebug = .true.
!
call myicopyso(iabij, r2, nioi3)
call myicopyso(fia, r1, nifvo)
!
maxvar = nifvo + nioi3 
iter = 0 
call printdnrm2('iabij', iabij, 3)
allocate(idoso(maxvar*10*isd), itmp(maxvar*isd))
!
!--iteration begins
!
do
iter = iter+1
write(6,*)
write(6,*)'=============in iteration:', iter
!
call cpu_time(tt(1)); call walltimeso(rt(1))
call denomt1(t1, r1, ie, 1)
call denomt2(t2, r2, ie, 1, 1, 1)
!
maxdim = 5
!
write(6,*)'begin dosodiis'
call cpu_time(tt1)
call dosodiis3(iter, t1, maxvar, maxdim, icore, &
              maxcore+maxvar, errdiis, idoso, itmp)
call cpu_time(tt2)
write(6,*)'time for dosodiis2',tt2-tt1
!
call effectt2(t2, t1,  icore, 1.d0, no, nv)
!
call calcene(t1, t2, fia, iabij, eccsd)
call calcene2(t1, t2, fia, iabij, eccsd2)
write(6,*)'CCSD energy', eccsd, eccsd2
eccsdio(iter) = eccsd2
write(6,*)'total energy', escf+eccsd, escf+eccsd2
ecc_i(iter) = eccsd2
!
call effectt2(t2, t1, icore, -1.d0, no, nv)
if(iter>29) exit
!
call cpu_time(tt(2)); call walltimeso(rt(2))
write(6,'(A, F10.3, F10.3)')'====time for energy     ', tt(2)-tt(1), rt(2)-rt(1)
!===========prapere for the begin 
call myicopyso(fab, fae, nifvv)
call myicopyso(fij, fmi, nifoo)
call myicopyso(fia, fme, nifvo)
call myicopyso(iijkl, wmnij, nioi1)
call myicopyso(fia, r1, nifvo)
call myicopyso(iabij, r2, nioi3)

call getallso_s(icore, 1, 'BIAJAAAS') 
call putallso_s(icore, 1, 'WMBEJAAA') 

call getallso_s(icore, 1, 'BIAJABAS') 
call putallso_s(icore, 1, 'WMBEJABA') 

call getallso_s(icore, 1, 'BIAJBBBS')
call putallso_s(icore, 1, 'WMBEJBBB') 

call getallso_s(icore, 1, 'BIAJBABS')
call putallso_s(icore, 1, 'WMBEJBAB') 

call getallso_s(icore, 1, 'BIAJAABS')
call putallso_s(icore, 1, 'WMBEJAAB') 

call itranspso(icore, icore(nioi46+1), & 
     vrta, popa, 1, vrtb, popb, 1, 1)
call putallso_s(icore(nioi46+1), 1, 'WMBEJBBA') 

call cpu_time(tt(3)); call walltimeso(rt(3))
write(6,'(A, F10.3, F10.3)')'====time for prepare    ', tt(3)-tt(2), rt(3)-rt(2)
!
!==========calc for the large terms 
!
!------term needed in r1 an r2 calculations 
call cpu_time(tt(4)); call walltimeso(rt(4))
call printdnrm2('fae', fae, 11)
call printdnrm2('ijkl', iijkl, 1)
call printdnrm2('ijka', iijka, 2)
call printdnrm2('abij', iabij, 3)
call printdnrm2('rivv', rivv, 111)
call printdnrm2('rivv', rivo, 222)
call printdnrm2('rivv', rioo, 333)

call calciabci_n5_ri_io(r1, r2, t1,t2, fae, icore, rivv, rivo, rioo)

call cpu_time(tt(5)); call walltimeso(rt(5))
write(6,'(A, F10.3, F10.3)')'====time for wmbej large', tt(5)-tt(4), rt(5)-rt(4)
call printdnrm2('fae', fae, 11)
call cpu_time(tt(5)); call walltimeso(rt(5))
!==============================
!-----terms which are independent 
!=============================
!
!==========abci======================
call cpu_time(tt(6)); call walltimeso(rt(6))
 call calciabci_ri(r1, r2, t1, t2, icore, rivv, rivo, rioo)
call cpu_time(tt(7)); call walltimeso(rt(7))
write(6,'(A, F10.3, F10.3)')'====time for iabci      ', tt(7)-tt(6), rt(7)-rt(6)
!
!==========abcd======================
call calciabcd_ri(r1, r2, t1, t2, icore, rivv)
call cpu_time(tt(8)); call walltimeso(rt(8))
write(6,'(A, F10.3, F10.3)')'====time for iabcd      ', tt(8)-tt(7), rt(8)-rt(7)
!
!=======    wmbej   ========================================
!
call cpu_time(tt(9)); call walltimeso(rt(9))
call calcwmbej_io(iijka, iabij, &
                t1, t2, r1, r2, &
                icore)
call cpu_time(tt(10)); call walltimeso(rt(10))
write(6,'(A, F10.3, F10.3)')'====time for wmbej      ', tt(10)-tt(9), rt(10)-rt(9)
!==============left============================
call calcleft_io(iijka, iabij, &
          t1, t2, r1, r2, &
          fae, fmi, fme, &
          wmnij, icore)
call cpu_time(tt(11)); call walltimeso(rt(11))
write(6,'(A, F10.3, F10.3)')'====time for left       ', tt(11)-tt(10), rt(11)-rt(10)
call printdnrm2('r1', r1, 22)
call printdnrm2('r2', r2, 3)
!if(igpuabcd==0) then 
!call myicopy_g(r1_c, icore, nifvo) 
!call axpyso(nifvo, 1.d0, icore, 1, r1, 1)
!call myicopy_g(r2_c, icore, nioi3) 
!call axpyso(nioi3, 1.d0, icore, 1, r2, 1)
!endif
!
!----
if(ldebug) then 
call printdnrm2('r1', r1, 22)
call printdnrm2('r2', r2, 3)
endif
!
write(6,'(A, F10.3, F10.3)')'====time for iteration  ', tt(11)-tt(1), rt(11)-rt(1)
iter_time(iter) = rt(11)-rt(1)
!
enddo
write(6,*)'the final energy is', Eccsd2
do i = 1, iter 
   write(6,'(A, I4, 2F18.12)')'ECCSD with ITER', i, Ecc_i(i), iter_time(i)
enddo
open(555,file='CCSD_ENERGY')
write(555,*) ESCF, Eccsd2
close(555)
!
      WRITE(6,99)
99    FORMAT(T3,'     Summary of iterative solution of CC equations ')
      WRITE(6,100)
100   FORMAT(6X,59('-'))
      WRITE(6,101)
101   FORMAT(T24,'Correlation',T46,'Total',/,T8,'Iteration',T26, &
             'Energy',T46,'Energy')
      WRITE(6,100)
      DO I=1,ITER
       WRITE(6,1000) I-1, Eccsdio(I), Eccsdio(I)+ESCF, ' DIIS '
1000   FORMAT(T10,I4,T19,F18.12,T39,F19.12,T59,A6)
      ENDDO
      WRITE(6,100)

      T1DIAG = 0.D0
      T1DIAG = nrm2so(irpdso(1,nva, noa, 1), T1, 1)**2
      T1DIAG = DSQRT(T1DIAG/dble(NOA))
      WRITE(6,*) '     T1 Diagnostic =', T1DIAG
!
      WRITE(6,103) ESCF
      WRITE(6,102) ECCSD2 + ESCF
103   FORMAT('        E(SCF)            = ',F20.12)
102   FORMAT('        E(CCSD)           = ',F20.12)

return
end 
